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The  approach  combining  image  analysis  techniques  and  artificial  neural  networks  is  proposed  here  for 
automatic  classification  of  mineral  inclusions  and  pores  in  archaeological  potsherds  using  optical  digital 
images.  Particularly,  the  automatic  identification  of  quartz,  calcareous  aggregates  and  secondary  porosity 
is  considered.  A  collection  of  both  plane  and  cross  polarised  light  images  acquired  via  a  digital  camera 
connected  to  optical  microscopy  in  transmitted  light  is  used.  Images  concern  Holocene  potsherds  (8900 
-4200  years  BP)  from  Takarkori  rock  shelter  archaeological  site  (SW  Libya,  Central  Sahara). 

The  adopted  methodology  involves  different  phases.  Firstly,  image  segmentation  is  carried  out  to 
isolate  regions  corresponding  to  the  interested  mineral  inclusions  and  pores.  A  segmentation  procedure 
based  on  mathematical  operators  is  customized  for  each  type  of  inclusions  and  for  pores.  Secondly, 
numerical  features  are  extracted  from  each  segmented  region,  thus  collecting  data  to  perform  automatic 
classification.  A  modular  classifier  is  considered  for  classification,  which  is  based  on  a  combination  of 
three  two-layer  feed-forward  neural  networks  that  are  trained  separately  to  recognise  each  class. 

Experimental  results  show  that  the  created  modular  classifier  provides  high  classification  accuracy  for 
both  inclusions  and  pores.  The  classifier  was  finally  applied  absent  the  image  analysis  phase  on  new 
samples  to  show  the  effectiveness  of  the  proposed  methodology. 

©  2014  Elsevier  Ltd.  All  rights  reserved. 


1.  Introduction 

In  ceramic  petrography,  optical  polarising  microscope  is 
routinely  used  to  study  archaeological  ceramics  via  observing  thin 
sections  of  these  materials.  In  this  way,  the  fabric— as  defined  by  the 
mineralogical  composition,  percentage,  shape,  size  distribution  and 
orientation  of  argillaceous  minerals,  non  plastic  inclusions  and 
pores— can  be  characterised  and  classified.  Thus,  the  definition  of 
raw  material  sources  and  selection,  forming  techniques,  chronol¬ 
ogy,  decorative  systems  and  firing  conditions  of  ancient  ceramics 
may  be  attained  (Maggetti,  1982;  Stoltman,  1989).  However, 
drawbacks  such  as  subjectivity  and  time-consuming  of  this  practice 
should  be  overcome.  Image  analysis  techniques  gave  an  important 
contribution  to  such  a  purpose.  Shape,  size,  roundness,  etc.  of 
minerals  or  grains  may  be  actually  automatically  measured  in  this 
way  and  used  to  characterise  rocks  (Hofmann  et  al.,  2011,  2013; 
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Holden  et  al.,  2009;  Perring  et  al.,  2004;  Smith  and  Beermann, 

2007)  as  well  as  archaeological  ceramics  and  mortars,  which  can 
be  actually  considered  as  artificial  rocks  (Caro  and  Di  Giulio,  2004; 
Casadio  et  al.,  2005;  Eramo  et  al.,  2014;  Livingood  and  Cordell, 
2009;  Middleton  et  al.,  1985;  Miriello  and  Crisci,  2006). 

Recently,  artificial  neural  networks  (ANNs)  coupled  to  image 
analysis  used  to  perform  automatic  classification  were  considered 
also  in  optical  mineralogy.  In  this  approach,  data  extracted  by  image 
analysis  relative  to  objects  of  interest  are  used  as  features  for  su¬ 
pervised  training  of  ANNs.  When  presented  with  an  unknown  ob¬ 
ject,  the  network  should  have  the  ability  to  assign  the  object  to  its 
class,  i.e.  to  finally  classify  that  object. 

Since  some  minerals  share  similar  properties  (particularly 
colour  and  brightness)  depending  on  the  interaction  with  light, 
image  analysis  based  on  optical  digital  images  generally  requires 
the  assessment  of  ad-hoc  hardware  and  processing  algorithms  to 
easily  distinguish  minerals  during  segmentation  (Fueten  and 
Goodchild,  2001;  Herwegh,  2000;  Soroushian  et  al.,  2003;  Wilson 
et  al.,  2007).  However,  although  other  image  acquisition  systems 
(e.g.  digital  X-ray  maps  by  scanning  electron  microscopy)  are  also 


A.  Aprile  et  al.  /  Journal  of  Archaeological  Science  50  (2014)  262—272 


263 


Fig.  1.  Example  of  thin  section  belonging  to  QC  fabric. 


used  (Knappett  et  al.,  2011;  Pret  et  al.,  2010),  image  analysis  based 
on  optical  microscopy  may  provide  a  wide  range  of  suitable  fea¬ 
tures  to  automatically  classify  minerals  and  textures.  Moreover, 
optical  microscopy  may  be  considered  to  be  low-cost  both  for 
sample  preparation  and  for  the  required  hardware,  especially  if 
combined  with  the  use  of  a  freeware  for  image  analysis. 

Despite  ANNs  coupled  to  optical  image  analysis  is  used  for 
automatic  classification  of  minerals  and  textures  of  rocks  (e.g. 
Baykan  and  Yilmaz,  2010;  Fueten  and  Mason,  2007;  Marmo  et  al., 
2005;  Singh  et  al.,  2010;  Thompson  et  al.,  2001),  to  our  knowl¬ 
edge  this  approach  has  not  been  applied  till  now  to  archaeological 
materials  such  as  pottery.  Actually,  ANNs  are  used  only  at  a  macro¬ 
scale  for  classification  of  potsherd  decorations  or  as  statistic 
method  in  chemiometric  studies  (Ma  et  al.,  2000;  Remola  et  al., 
1996;  Stanco  et  al.,  2011).  However,  as  well  as  for  rocks,  the 
learning  and  generalization  abilities  of  neural  classifiers  may  be 
useful  to  “translate”  heterogeneity  of  archaeological  materials  to 
attain  their  classification  and  characterisation. 

Such  an  approach  is  thus  proposed  in  this  paper  for  archaeo¬ 
logical  pottery,  specifically  to  face  the  problem  of  automatic  iden¬ 
tification  of  mineral  inclusions  and  porosity.  Particularly,  quartz, 
calcareous  aggregates  and  secondary  porosity  are  considered. 

The  adopted  methodology  involved  different  phases.  Digital 
plane  and  cross  polarised  light  images  were  both  acquired  by  thin 
sections.  Image  processing  was  then  performed  to  segment  regions 
corresponding  to  the  interested  inclusions  and  pores.  Statistical  and 
region  features  were  hence  extracted  by  each  segmented  region  to 
be  used  for  neural  network  training.  Particularly,  a  modular  archi¬ 
tecture  is  used  based  on  a  combination  of  three  neural  networks 
each  trained  to  recognise  one  class.  The  created  modular  neural 
classifier  was  finally  applied  absent  the  image  analysis  phase  on 
new  samples  to  evaluate  its  practical  application  for  these  types  of 
mineral  inclusions  and  for  pores. 

2.  Materials 

Some  Holocene  potsherds  (8900-4200  BP)  belonging  to  the 
Takarkori  archaeological  site  (SW  Libya,  Central  Sahara)  are 


considered  here  as  case  study.  A  complete  description  of  the 
archaeological  context  can  be  found  in  Eramo  et  al.  (2014). 

Petrographic  observations  distinguished  three  fabrics  (named  Q 
QC  and  QF)  for  these  potsherds  according  to  different  prevailing 
clastic  constituents.  The  QC  fabric  is  characterised  by  calcareous 
aggregates  and  quartz  grains  as  mineral  inclusions.  The  matrix 
presents  high  optical  density  due  to  carbonaceous  matter,  in  some 
cases  partially  oxidised  near  the  external  surface  and  secondary 
pores  related  to  the  partial  combustion  of  vegetal  inclusions  and 
organic  matter,  as  well  as  to  sintering  occurred  during  firing.  The  QF 
fabric  is  rather  characterised  by  potassium  and  plagioclase  feld¬ 
spars  along  with  angular  or  sub-angular  grains  of  quartz  as  mineral 
inclusions.  The  Q  fabric  has  a  very  similar  petrofacies  compared  to 
QC  because  characterised  only  by  quartz  inclusions.  These  two 
fabrics  were  already  studied  by  using  image  analysis.  Since  of  their 
mineralogical  similarity,  grain-size  and  morphometric  data  ob¬ 
tained  in  this  way  for  quartz  and  calcareous  aggregate  inclusions 
were  useful  to  determine  the  common  fluviatile  origin.  Moreover,  a 
compatibility  with  wadi  and  swamp  sedimentary  environments 
present  in  the  area  was  found,  which  confirmed  the  local  prove¬ 
nance  of  the  raw  materials  (Eramo  et  al.,  2014). 

The  14  thin  sections  of  potsherds  belonging  to  fabric  QC  were 
used  to  assess  the  proposed  methodology  (Fig.  1). 

3.  Image  analysis 

The  image  analysis  phase  is  devoted  to  create  a  set  of  labelled 
data  to  be  used  for  the  supervised  learning  of  the  ANNs.  Particu¬ 
larly,  in  this  phase,  the  knowledge  about  the  interested  classes, 
namely  quartz,  calcareous  aggregates  and  pores,  was  used  to 
segment  images  into  the  corresponding  regions.  After  segmenta¬ 
tion,  feature  extraction  is  performed  so  as  to  describe  each  region 
by  a  vector  of  numerical  features.  Precisely,  a  set  of  labelled  samples 
(i.e.  the  segmented  regions)  for  each  class  c  =  1  ...  I<  (in  our  case 
I<  =  3)  has  to  be  created,  each  described  by  a  feature-vector  to  be 
used  for  training  the  neural  classifier. 

To  perform  image  analysis,  four  pairs  of  RGB  digital  images  each 
composed  by  a  plane  and  a  cross  polarised  light  image  (denoted  in 
the  following  respectively  as  P  and  XP)  of  1920  x  2560  pixel  size 
were  acquired  per  thin  section  at  a  magnification  of  2.5  x  (1.315  pm / 
pixel)  and  then  saved  in  TIFF  format  (Fig.  2).  A  ZEISS  Axioskop  40  Pol 
petrographic  microscope  equipped  with  a  Nikon  DS-5MC  digital 
camera  with  2/3”  CCD  (Peltier  cooling  system)  was  used  for  acqui¬ 
sition.  The  four  pairs  of  P  and  XP  images  were  acquired  by  four 
different  locations  of  each  thin  section,  defined  by  maintaining  a 
distance  of  1  cm  between  the  centre  of  each  frame  by  the  successive 
one  and  moving  by  one  side  to  the  other  of  the  thin  section. 


Fig.  2.  Plane  and  cross  polarised  light  digital  images  acquired  by  optical  microscopy  (example  shown  refers  to  sample  named  tl<7). 
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Fig.  3.  Segmentation  of  quartz. 


3  A.  Segmentation 

To  segment  regions  corresponding  to  the  three  classes,  both  the 
acquired  P  and  XP  images  were  considered.  The  freeware  open- 


source  image  analysis  software  Image /  (version  1.42  and  later; 
Rasband,  1997-2012)  was  used  to  such  a  purpose. 

The  manipulation  of  both  kinds  of  images  guarantees  a  better 
identification  of  minerals  and  the  detection  of  characteristics  useful 
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for  segmentation.  Particularly,  as  highlighted  by  Whitbread  (1991 ), 
inclusion  morphology  is  more  clearly  shown  by  the  P  images,  while 
their  mineralogical  distinction  is  more  easily  achieved  considering 
the  interference  colours  revealed  by  the  XP  ones.  The  segmentation 
phase  was  primarily  based  on  the  use  of  mathematical  operators  to 
emphasize  useful  characteristics  of  the  interested  inclusions  and 
pores.  Then,  the  Isodata  thresholding  method  (Ridler  and  Calvard, 
1978)  was  finally  applied.  In  this  way,  a  customized  segmentation 
procedure  was  developed  for  each  class  in  order  to  derive  corre¬ 
sponding  targeted  samples. 

3.2.2.  Segmentation  of  quartz 

For  segmentation  of  quartz,  both  P  and  XP  images  were 
manipulated  together.  Actually,  some  quartz  grains  are  not  visible 
using  only  XP  images.  At  the  same  time,  quartz  grains  are  not  easily 
distinguishable  from  pores  if  only  P  images  are  considered,  since  in 
these  images  the  colours  of  pores  (in  the  range  of  white  values)  are 
similar  to  the  colours  of  quartz. 


Particularly,  P  and  XP  images  were  firstly  combined  using  the 
addition  operator  to  obtain  a  new  image  where  all  the  regions  cor¬ 
responding  to  quartz  are  morphologically  well  distinguishable  and 
pixel  brightness  was  minimized  to  have  a  better  edge  definition 
(Fig.  3A).  However,  the  obtained  image  still  showed  the  same  range  of 
colours  for  both  quartz  and  pores.  To  avoid  this  effect,  the  XP  image 
was  subtracted  by  the  P  one  using  the  subtraction  operator  (Fig.  3B). 
This  procedure  was  inspired  by  the  method  proposed  by  Marinoni 
et  al.  (2005)  for  an  experimental  mortar  thin  section  study  via  im¬ 
age  analysis  aimed  to  isolate  porosity.  Using  this  procedure,  many 
quartz  regions  were  distinguishable  from  pores  since  they  were 
characterized  by  darker  colours  (grey-brownish  areas  of  the  image). 

Lastly,  to  remove  quartz  regions  that  were  still  confused  with 
pores,  the  exponential  mathematical  operator  was  applied  to  the 
previously  obtained  image,  emphasizing  high  pixel  values  (Young 
et  al.,  1998)  corresponding  to  pores  (Fig.  3C).  In  this  way,  segmen¬ 
tation  of  quartz  grains  was  obtained  by  progressively  emphasizing 
the  original  range  of  colour  values  of  pores. 


Fig.  4.  Segmentation  of  calcareous  aggregates. 


266 


A.  Aprile  et  al.  /  Journal  of  Archaeological  Science  50  (2014)  262—272 


Successively,  the  obtained  image  was  subtracted  by  the  image 
resulting  by  the  addition  operation  in  the  first  step  in  order  to 
finally  isolate  quartz  grains  (Fig.  3D).  Then,  it  was  converted  to  8-bit 
greyscale  format  (Fig.  3E)  to  be  finally  thresholded  applying  the 
Isodata  method  (Fig.  3F).  The  resulting  binary  map  was  then  pro¬ 
cessed  using  the  pixel  eraser  and  the  hole  filling  morphological 
operators  to  obtain  the  final  binary  image  containing  quartz  regions 
(Fig.  3G). 

3.1.2.  Segmentation  of  calcareous  aggregates 

To  segment  calcareous  aggregates,  the  differences  in  colours  and 
texture  with  respect  to  quartz  were  considered.  Actually,  the 
calcareous  aggregates  show  a  sponge-like  texture  with  less  sharpen 
borders  than  quartz.  Also,  they  show  a  colour  range  characterized 
by  a  prevalence  of  R  (red)  and  G  (green)  combination  (i.e.  a 
particular  light-to-dark  brown  colour).  Moreover,  a  lower  contrast 
was  observed  with  respect  to  quartz.  These  differences  are  more 
clearly  visible  in  XP  images,  therefore  segmentation  was  performed 
using  these  images. 


Firstly,  the  quartz  binary  image  previously  obtained  (Fig.  3G) 
was  inverted  (i.e.  the  white  colour  was  assigned  to  the  inclusions) 
and  filtered  using  the  Gaussian  Blur ,  the  Median  and  the  Maximum 
filters  (Fig.  4A).  These  filters  were  used  in  order  to  guarantee  a 
better  overlap  of  the  edges  of  grains  in  the  binary  image  with 
respect  to  the  corresponding  ones  in  XP  image  and  selected  ac¬ 
cording  to  characteristics  of  each  processed  image.  The  filters  were 
parameterized  according  to  the  default  setting  of  the  software, 
namely  sigma  =  2.0  for  the  Gauss  filter  and  radius  =  2.0  for  the 
Median  and  the  Maximum  filters. 

Next,  the  CLAHE  (Contrast  Limited  Adaptive  Fiistogram  Equal¬ 
ization)  method  (Zuiderveld,  1994)  was  used  for  contrast 
enhancement  of  XP  image  to  better  distinguish  calcareous  aggre¬ 
gates  from  quartz  (Fig.  4B).  The  quartz  binary  image  initially 
inverted  and  filtered  was  then  subtracted  by  corresponding  XP 
image  previously  enhanced  (Fig.  4C).  The  resulting  image  was 
finally  converted  in  the  8-bit  greyscale  format  (Fig.  4D)  and  then 
thresholded  to  produce  the  binary  image  of  calcareous  aggregates 
(Fig.  4E). 


Fig.  5.  Segmentation  of  pores. 
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3.1.3.  Segmentation  of  pores 

For  the  segmentation  of  pores,  only  P  images  were  considered 
because  they  more  clearly  show  morphology  and  colours  of  pores. 
Firstly,  both  the  quartz  and  the  calcareous  aggregate  binary  images 
(Figs.  3G  and  4E)  were  inverted  and  filtered  using  the  same  filters 
previously  cited  to  guarantee  a  better  overlap  with  the  P  image  in 
the  next  steps  (Fig.  5 A).  Then,  each  image  was  subtracted  by  the  P 
image.  The  obtained  image  (Fig.  5B)  was  thus  converted  in  the  8-bit 
greyscale  format  (Fig.  5C).  The  binary  image  showing  pore  regions 
was  produced  after  thresholding  and  application  of  the  median 
filter  Despeckle  (Fig.  5D)  implemented  in  Image j  to  reduce  noise, 
which  exactly  replaces  each  pixel  value  with  the  fifth  median  value 
in  its  3  x  3  neighbourhood  (Ferreira  and  Rasband,  2010-2012). 

3.2.  Feature  extraction 

To  evaluate  characteristics  of  inclusions  and  pores  (i.e.  colours, 
shape,  orientation,  etc.)  and  thus  to  accomplish  the  classification 
task,  the  neural  classifier  requires  numerical  features  describing 
quantitatively  these  characteristics.  As  in  other  similar  works 
(Baykan  and  Yilmaz,  2010;  Singh  et  al.,  2010),  statistical  and  region 
features  were  considered  here.  Such  features  were  computed  by 
using  the  Image ]  plug-in  Particles8_Plus  (Landini,  2008). 

Statistical  features  are  descriptors  of  intensity  histogram  of  the 
image  or  of  the  considered  segmented  regions  related  to  moments 
of  distribution  and  standard  statistics  of  colours  or  grey  tones 
(Gonzalez  and  Woods,  2002;  Press  et  al.,  1992).  About  region  fea¬ 
tures,  these  parameters  may  be  defined  considering  the  external 
descriptors  proposed  by  Gonzalez  and  Woods  (2002)  for  regions  in 
an  image.  Flowever,  different  denominations  are  reported  here 
according  to  the  ones  used  by  the  adopted  plug-in. 

We  extracted  features  by  each  singular  region  of  quartz, 
calcareous  aggregate  and  pore  as  obtained  by  segmentation,  which 
refer  to  the  corresponding  areas  of  P  and  XP  images.  Particularly, 
the  statistical  features  were  computed  for  each  R  (red),  G  (green) 
and  B  (blue)  colour  channel  related  to  the  distribution  histogram  of 
the  corresponding  images.  A  39  feature-vector  for  each  region  was 
thus  provided  to  perform  automatic  classification.  The  extracted 
features  are  reported  and  briefly  described  in  able  1. 

4.  The  modular  neural  classifier 

To  distinguish  among  quartz,  calcareous  aggregates  and  pores,  a 
modular  architecture  where  each  classifier  is  trained  to  classify  one 
of  the  interested  classes  (Jacobs  and  Jordan,  1993;  Sarkar,  2000)  was 

adopted  here.  In  this  architecture,  the  classification  problem  is 
decomposed  into  a  number  of  sub-problems  addressed  by  modules 
that  operate  simultaneously  on  the  same  inputs  without  commu¬ 
nicating  each  other.  A  combining  scheme  recombines  the  outputs 
of  the  modules  to  provide  the  final  classification.  The  use  of  a 
modular  classifier  can  help  to  improve  the  classification  accuracy 
having  the  advantages  of  quick  convergence,  parallel  training  and 
better  treatment  of  uncertainties  (Briem  et  al.,  2002;  Sarkar,  2000). 

Although  different  classification  methods  can  be  used  to  solve 
each  sub-problem,  the  use  of  ANNs  for  each  module  was  preferred 
here.  ANNs  are  computational  modelling  tools  that  emulate  the 
biological  (human  or  animal)  central  nervous  system  in  performing 
a  particular  task.  An  ANN  is  structured  as  a  densely  interconnected 
system  of  “artificial  neurons”,  i.e.  simple  processing  elements, 
forming  input,  hidden  and  output  layers  capable  of  acquisition, 
parallel  data  computation  and  knowledge  representation  (Bishop, 
1996).  ANNs  are  particularly  suited  in  pattern  classification  or  ob¬ 
ject  recognition  problems  and  when  no  algorithm  or  model  is 
available  to  solve  a  problem;  rather  a  large  amount  of  data  are 
available  as  examples  to  be  used  for  supervised  learning. 


Table  1 

Features  used  in  this  study. 


Statistical 

features 

Definition 

Region  features 

Definition 

Mode 

Most  occurring  pixel  value 

Solidity 

Area/convex  area 

Median 

Median  pixel  value 

Concavity 

Convex  hull-area 

Average 

Mean  pixel  value 

Rectangularity 

Area/area 
bounding  box 

Variance 

Mean  square  deviation 
of  pixel  values  from  average 

Standard 

Standard  deviation 

deviation 

of  pixel  values 

Skewness 

Degree  of  asymmetry 
of  pixel  value  distribution 

Kurtosis 

Peakedness  of  pixel 
value  distribution 

Average 

Spread  of  pixel  values 

deviation 

from  average 

Integrated 

Sum  of  pixel  values 

density 

Min 

Minimum  pixel  value 

Max 

Maximum  pixel  value 

Entropy 

Degree  of  variability 
of  pixel  values 

In  the  following,  the  neural  network  model  used  to  implement 
the  neural  modules  and  the  scheme  for  combining  their  outputs  are 
described. 


4.1.  Creation  of  modules 

As  diffusely  used  for  pattern  recognition  problems  (Arel,  2012; 
Carson  et  al.,  1995;  Giiler  and  Ubeyli,  2005;  Yildirim  et  al.,  2004), 

a  Multilayer  Perceptron  Neural  Network  (MLPNN)  model  was 
adopted  here  to  implement  each  module.  It  is  a  feed-forward 
network  generally  consisting  of  an  input  layer  with  neurons  rep¬ 
resenting  input  variables  describing  the  problem,  an  output  layer 
whose  neurons  correspond  to  variables  describing  what  is  being 
learnt  and  one  or  more  hidden  layers  having  neurons  capable  to 
capture  the  non  linearity  of  the  data.  The  set  of  outputs  signals  in 
the  output  layer  constitutes  the  response  of  the  network  to  the 
pattern  provided  as  input  (Bishop,  1996;  Haykin,  1999). 

Particularly,  in  this  work  a  two-layer  MLPNN  with  39  inputs  was 
used.  One  output  is  expected  in  the  output  layer  since  each  module 
is  created  to  classify  one  of  the  considered  classes  (Fig.  6). 

The  number  of  neurons  in  the  hidden  layer  was  ranged  from  5  to 
80  using  a  rate  of  5  neurons.  In  this  way,  a  total  amount  of  16 
different  two-layer  MLPNN  topologies  were  trained  for  each  class. 
The  MatLab®  NPRTool  (Neural  Network  Pattern  Recognition  Tool) 
v.7.10.0  (R2010a)  was  used  to  such  a  purpose. 

To  create  the  modules,  the  original  three-class  classification 
problem  was  primarily  divided  into  three  two-class  problems 
following  the  one-vs-rest  task  decomposition  scheme  (Anand  et  al., 
1995;  Chen  and  You,  1993).  According  to  this  scheme,  each  neural 
module  is  used  to  discriminate  patterns  belonging  to  one  class  from 
patterns  belonging  to  the  remaining  classes.  Therefore,  the  i-th 
module  or  classifier  is  trained  by  using  all  patterns  of  class  c  as 
positive  samples  T+  and  the  remaining  patterns  as  negative  samples 
Tf.  The  training  set  for  each  of  the  two-class  problems  was  theo¬ 
retically  defined  by 

T  c  =  {xi,yic}^1  for  c  =  1...K 

where  Nc  is  the  number  of  training  data,  I<  is  the  number  of  classes 
(here  I<  =  3),  X;  is  the  input  vector  of  the  i-th  training  sample  and  yIC 
is  the  desired  output  for  the  i-th  training  sample  which  is  defined  as 
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Input-vector  xi 


Input-vector  X2 


Input-vector  X3 


Input-vector  X39 


►  output 


INPUT  HIDDEN  OUTPUT 

LAYER  LAYER  LAYER 

Fig.  6.  The  two-layer  MLPNN  network  used  in  this  study. 


_  f  1  if  Xj  belongs  to  class  c 
^lc  ~  \  0  otherwise 

Particularly,  the  subset  of  positive  samples  for  class  c  is 

T+  =  {xl-,yIC}^1  where  N £  is  the  number  of  samples  representative 
of  class  c,  Xi  is  the  feature-vector  of  the  i-th  sample  and  yIC  is  the 
desired  output  for  the  i-th  sample  which  is  defined  as  yIC  =  1.  The 
subset  of  negative  samples  for  class  c,  uniformly  taken  from  the  sets 

TjJ-,  with  k^=c,  is  denoted  as  T^  =  {Xj,ylc}J^1  where  N k  is  the 
number  of  samples  that  do  not  represent  class  c,  Xj  is  yet  the 
feature-vector  of  the  i-th  sample  and  ylc  is  the  desired  output  for 
the  i-th  sample  which  is  defined  as  yIC  =  0.  Precisely,  to  create  a 
neural  module  capable  of  recognising  each  class  c,  the  creation  of  a 

dataset  Tc  =  T+uT~  =  {xz-,ylc}-l1  is  needed,  where  N  =  N+  +  Nk  is 
the  total  number  of  samples. 

Given  the  training  set  Tc,  a  binary  classifier  has  to  be  created 
which  is  capable  to  classify  each  input  pattern  in  the  correct  class, 
i.e.  to  discriminate  one  class  from  the  remaining  two  classes. 

To  perform  training,  each  dataset  was  randomly  split  here  into  a 
70%  training  set,  a  15%  test  set  and  a  15%  validation  set.  A  supervised 
learning  scheme  called  scaled  conjugated  backpropagation  algo¬ 
rithm  was  adopted  (Moller,  1993).  Training  was  stopped  when  the 
minimum  error  on  the  validation  set  was  reached  or  when  a 
maximum  of  1000  epochs  was  achieved.  To  achieve  more  reliable 
results,  10  training  runs  were  performed  for  each  of  the  16  network 
topologies  starting  from  a  different  initial  configuration  of  weights. 

The  classification  error  on  the  test  set,  defined  as  the  probability 
of  error  in  classifying  new  objects  (Wanas  et  al.,  1999),  was 
considered  to  evaluate  the  classification  accuracy  of  the  trained 
networks.  Precisely,  the  classification  accuracy  is  calculated  as  the 
number  of  correctly  classified  samples  over  the  total  number  of  the 
considered  ones,  that  is  the  sum  of  the  true  positive  samples  and 
the  true  negative  samples  over  the  total  number  of  the  considered 
ones  (Fawcett,  2006). 

To  select  the  optimal  neural  network  topology  to  be  used  as 
neural  module  in  the  modular  architecture  for  each  class,  two 
different  criteria  were  investigated.  According  to  the  former,  the 
network  topology  providing  the  lowest  classification  error  on  the 
test  set  among  the  160  trained  networks  was  selected.  According 


1  The  training  error  generally  decreases  as  training  progresses  while  the  valida¬ 
tion  set  error  only  decreases  up  to  a  certain  point  after  which  it  increases  and  this 
generally  happens  when  an  over-fitting  of  the  data  is  reached  by  the  network. 


to  the  latter,  the  network  topology  providing  the  lower  mean 
value  of  the  classification  error  on  the  test  set  computed  on  the 
corresponding  10  training  runs  was  primarily  considered.  Then, 
the  configuration  which  returned  the  lower  classification  error 
among  the  10  training  run  recorded  for  that  topology  was  finally 
chosen. 


4.2.  Recombination  of  modules 

After  each  of  the  two-class  problems  was  learnt  by  a  MLPNN 
module,  an  output  /c(x*)  of  the  c-th  classifier  for  input  pattern  Xj  is 
returned  and  a  function  fc  is  created  so  as  to  discriminate  the  pat¬ 
terns  of  class  c  =  1  ...  I<  (with  I<  =  3)  from  those  of  the  other  classes. 
We  created  the  modular  classifier  by  applying  simultaneously  the 
three  neural  modules  to  an  input  pattern  or  assessment  dataset  and 
by  combining  their  outputs  in  order  to  derive  the  final 
classification. 

Several  combination  schemes  for  modular  classifiers  have  been 
proposed  in  the  literature  (Lu  and  Ito,  1999;  Sarkar,  2000).  Partic¬ 
ularly,  we  applied  a  standard  recombination  scheme  based  on  the 
winner-take-all  strategy,  that  works  as  follows.  We  denote  by  /(x) 
the  actual  output  vector  of  the  modular  classifier  for  the  three-class 
classification  problem.  It  can  be  written  as 

/(X)  =  [fi(x),/2(x),/3(x)]T 

Since  the  outputs  /c(x/)  returned  by  each  classifier  can  be  also 
interpreted  as  confidence  values,  a  possible  output  value  of  0.9,  for 
example,  means  that  a  pattern  is  assigned  to  the  positive  class  with 
a  high  confidence  while  a  possible  value  of  0.3,  on  the  other  hand, 
means  that  the  pattern  should  be  considered  negative  with  reser¬ 
vations.  Namely,  the  larger  the  value,  the  more  confidently  the 
point  Xi  belongs  to  the  class  c.  Thus,  given  an  unknown  input 
pattern  x*,  the  combining  scheme  takes  as  final  class  the  class  with 
the  highest  confidence  value.  Precisely,  the  modular  classifier  as¬ 
signs  x*  to  class  c  if  the  following  condition  holds: 

|/c(x*)  —  1|  <  5  and  |//<(x*)|  <5  for  k^c  (1) 

where  <5  is  a  real  number,  which  denotes  the  error  tolerance  (in  this 
work  5  is  set  to  0.5). 

For  example, /(x*)  =  [0.9,  0.4,  0.1  ]  means  that  pattern  x*  belongs 
to  quartz, /(x*)  =  [0.3,  0.7,  0.2]  means  that  pattern  x*  belongs  to 
calcareous  aggregates  and /(x*)  =  [0.3,  0.3,  0.8]  means  that  pattern 
x*  belongs  to  pores  (Fig.  7). 
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Input  pattern 
x* 


COMBINATION 
SCHEME 


Kx-)=[ Mx),/>(x)./*(x)]r 


/(x*)=[l,0,0]=class  1 
/tx*)=[0.1,0]=class2 
/(x’)=[0,0,l]=class  3 


Fig.  7.  Modular  architecture  adopted  for  classification. 
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Feature 
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Modular  neural 
classifier 


I  *)  class 


Fig.  8.  Application  of  the  created  modular  neural  classifier  absent  the  image  analysis  phase. 


When  condition  (1)  is  not  satisfied,  the  classification  of  x*  is 
considered  ambiguous.  For  example  /(x*)  =  [0.8,  0.9,  0.2]  means 
that  x*  is  recognised  as  belonging  to  both  quartz  and  calcareous 
aggregates.  In  this  case,  the  class  of  x*  is  decided  based  on  the 
winner-take-all  strategy  namely  the  pattern  x*  is  assigned  to  class  c 
where /c(x*)  =  max  /fc(x*),  i.e.  the  class  with  the  highest  response. 

k=t...I< 

Once  properly  assessed,  the  created  modular  classifier  can  be 
used  to  recognise  the  class  of  new  samples  absent  the  image 
analysis  that  did  the  initial  classification.  In  this  scenario,  the  user 
simply  needs  to  manually  select  a  ROI  (i.e.,  region  of  interest)  cor¬ 
responding  to  any  inclusion  or  pore  he  wants  to  know  the  class  in 
an  image  of  a  new  thin  section  of  any  potsherd  and  run  the  software 
that  computes  the  39  numerical  features  starting  from  that  ROI.  The 
corresponding  feature-vector  will  be  then  used  as  input  to  the 
modular  classifier,  which  will  finally  provide  the  class  of  that  in¬ 
clusion  or  pore.  Fig.  8  illustrates  this  process. 

5.  Results 

In  this  work,  a  total  of  41,085  regions  corresponding  to  quartz, 
9313  regions  corresponding  to  calcareous  aggregates  and  35,258 
regions  corresponding  to  pores  were  finally  obtained  by  the 
adopted  segmentation  procedure  (Section  3). 

Then,  in  order  to  avoid  the  problem  of  unbalanced  data  for  the 
three  classes,  the  number  of  samples  to  be  used  in  the  training  set 
for  each  class  was  established  by  considering  the  number  of 
available  samples  for  the  less  abundant  class  (i.e.  the  calcareous 
aggregates).  Thus,  half  of  these  samples  as  well  as  an  equal  number 
for  the  remaining  classes  were  considered  for  training. 


2  Recall  that  the  term  sample  refers  to  each  segmented  region  with  its  corre¬ 
sponding  39  feature-vector. 


The  remaining  samples  for  each  class  were  used  to  create  an 
assessment  dataset  to  test  the  modular  classifier,  which  was  finally 
composed  by  3584  unknown  samples,  including  respectively  1339 
samples  of  quartz,  1089  of  calcareous  aggregates  and  1156  of 
pores. 

Classification  errors  on  the  test  set  provided  by  the  trained 
networks  for  each  of  the  three  classes  are  shown  in  Fig.  9. 

According  to  the  first  selection  criterion  (Section  4.1 ),  the  lowest 
error  for  class  of  quartz  was  achieved  by  a  39  x  25  x  1  network,  for 
calcareous  aggregates  by  a  39  x  5  x  1  network  and  for  pores  by  a 
39  x  5  x  1  network  (  Fable  2).  These  selected  neural  classifiers, 
combined  into  the  modular  classifier,  finally  achieved  a  classifica¬ 
tion  accuracy  of  90.1%  when  tested  on  the  assessment  dataset. 

On  the  base  of  the  second  criterion,  the  lowest  mean  error  was 
recorded  for  the  39  x  45  x  1,  the  39  x  20  x  1  and  the  39  x  40  x  1 
network  topologies  respectively  for  class  of  quartz,  calcareous  ag¬ 
gregates  and  pores,  as  shown  in  Table  3.  The  three  selected  classi¬ 
fiers  were  combined  then  to  form  a  second  modular  classifier,  that 
provided  a  classification  accuracy  of  89.8%  on  the  same  assessment 
dataset. 

Since  the  two  created  modular  classifiers  provided  similar 
classification  accuracy,  the  accuracy-simplicity  trade-off  criterion 
was  considered  to  finally  select  the  best  one.  That  is,  the  first 
modular  classifier  was  actually  chosen  because  it  turns  out  to  have 
a  simpler  topology  (i.e.,  simpler  network  modules)  than  the  second 
one. 

As  shown  in  the  corresponding  confusion  matrix  (Fig.  10), 
classification  accuracy  of  the  final  modular  classifier  was  higher  for 
calcareous  aggregates  and  quartz,  while  it  was  lower  for  pores.  The 
quartz  was  mainly  confused  with  the  calcareous  aggregates  (70 
errors  over  76)  as  well  as  the  calcareous  aggregates  and  the  pores 
were  confused  with  the  quartz  (respectively  36  errors  over  45  and 
178  errors  over  235). 
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Fig.  9.  Classification  error  on  test  set  provided  by  the  trained  neural  networks. 


Table  2 

Neural  classifiers  selected  according  to  the  first  criterion. 


Neural  module 

Topology 

Classification  error 

Epochs 

Quartz 

39  x  25  x  1 

2.41 

76 

Calcareous  aggregates 

39  x  5  x  1 

0.58 

58 

Pore 

39  x  5  x  1 

0 

77 

Finally,  Fable  4  shows  the  results  of  classification  as  performed  by 
the  final  modular  classifier  absent  the  image  analysis  phase  on  9  new 
regions  of  inclusions  and  pores  randomly  selected  by  new  images  of 
potsherds.  As  confirmed  by  the  expert,  the  modular  classifier 
correctly  classified  all  these  new  regions  into  the  correct  class. 

Moreover,  the  final  modular  classifier  was  compared  with  a 
three-class  monolithic  classifier  (namely  capable  to  contemporarily 


Table  3 

Neural  classifiers  selected  according  to  the  second  criterion. 


Neural  module 

Topology 

Classification  error 

Epochs 

Quartz 

39  x  45  x  1 

3.29 

85 

Calcareous  aggregates 

39  x  20  x  1 

1.45 

53 

Pore 

39  x  40  x  1 

0 

31 

1263 

35.3% 

36 

1.0% 

178 

5.0% 

85.5% 

14.5% 

70 

2.0% 

1044 

29.1% 

57 

1.6% 

89.2% 

10.8% 

6 

0.2% 

9 

0.2% 

921 

25.7% 

98.4% 

1.6% 

94.3% 

5.7% 

95.7% 

4.3% 

79.7% 

20.3% 

90.1% 

9.9% 

1  2  3 


Target  class 


Fig.  10.  Confusion  matrix  of  modular  classifier.  A  classification  accuracy  of  90.1%  is 
reached  by  modular  classifier  on  3584  unknown  samples,  which  is  higher  for  calcar¬ 
eous  aggregates  and  quartz  (respectively  94.3%  and  95.7%),  while  it  is  lower  for  pores 
(79.7%).  Quartz  is  mainly  confused  with  calcareous  aggregates  (70  errors  over  76)  as 
well  as  calcareous  aggregates  and  pores  are  confused  with  quartz  (respectively  36 
errors  over  45  and  178  errors  over  235). 
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Table  4 

Results  of  classification  performed  by  the  created  modular  classifier  absent  the  image  analysis  phase.  As  confirmed  by  the  expert,  the  modular  classifier  correctly  classified  all 
the  new  considered  regions  of  inclusions  and  pores  (ROI),  randomly  selected  by  new  images  of  potsherds,  into  the  correct  class. 

ROI  Class  ROI  Class  ROI  Class 


distinguish  among  quartz,  calcareous  aggregates  and  pores)  prop¬ 
erly  selected  from  several  two-layer  network  topologies  with  39 
inputs  and  3  outputs  created  using  the  same  training  conditions, 
which  finally  had  a  39  x  10  x  3  topology.  As  it  can  be  seen  in 
Table  5,  the  modular  classifier  reached  the  better  classification 
performance  towards  the  monolithic  network  when  tested  on  the 
same  assessment  dataset,  except  for  quartz.  Moreover,  for  both 
classifiers,  the  higher  number  of  errors  correspond  to  the  class  of 
pores,  which  has  the  lowest  classification  accuracy.  The  lower 
performance  of  the  monolithic  classifier  is  also  shown  by  the 
number  of  epochs  after  which  it  converged  to  a  solution  (namely 
106  epochs)  which  is  higher  than  the  ones  of  the  classifiers 
composing  the  modular  classifier  shown  in  Table  2. 

6.  Discussion  and  conclusions 

In  this  paper,  a  methodology  for  automatic  classification  of  two 
types  of  mineral  inclusions  (i.e.,  quartz  and  calcareous  aggregates) 
and  of  pores  in  archaeological  potsherds  is  presented  which  com¬ 
bines  image  analysis  and  ANNs.  Segmentation  of  inclusions  and 
pores  was  performed  using  mathematical  operators  on  both  P  and 
XP  images.  An  automatic  thresholding  method  was  finally  coupled 
to  improve  region  detection.  Then,  rather  than  using  a  monolithic 
classifier,  a  modular  architecture  that  combines  as  many  neural 
classifiers  as  the  number  of  classes  was  adopted  for  automatic 
classification. 

At  present,  our  preliminary  results  are  encouraging  and  indicate 
that  reliable  features  obtained  by  taking  into  account  optical 
characteristics  of  the  interested  classes  and  a  specific  neural 
network  architecture  can  provide  a  good  automatic  classification. 
Particularly,  the  better  classification  accuracy  achieved  for  calcar¬ 
eous  aggregates  can  be  referred  to  their  specific  colours  and  texture 
(see  Section  3.1.2)  which  probably  led  to  the  extraction  of  more 
distinctive  features  for  this  class.  The  confusion  between  quartz  and 


calcareous  aggregates  may  be  referred  to  some  samples  of  calcar¬ 
eous  aggregates  which  appeared  more  light-blue  coloured  with 
respect  to  their  common  brownish  colour,  thus  resembling  typical 
quartz  colouring.  Meanwhile,  the  confusion  between  quartz  and 
pores  may  be  related  to  some  samples  of  quartz  which  appeared 
lighter  assuming  similar  colours  as  the  pores.  However,  when 
tested  absent  the  image  analysis  phase,  the  modular  classifier 
assigned  all  the  new  samples  to  the  correct  class. 

Moreover,  the  simulation  results  show  that  the  classification 
accuracy  of  the  proposed  modular  classifier  is  about  the  same  as  the 
corresponding  monolithic  one.  However,  while  the  classification 
performance  of  the  modular  classifier  may  be  improved  choosing 
other  combinations  of  classifiers,  the  monolithic  classifier  topology 
remains  fixed.  This  means  that  another  classifier  for  the  class  of 
pore  may  be  selected  (also  adopting  other  selection  criteria)  and 
combined  into  a  new  architecture,  trying  to  reach  better  classifi¬ 
cation  performance.  Actually,  the  use  of  a  modular  classifier  pro¬ 
vides  the  possibility  to  grow  up  the  original  architecture  with 
additional  modules  as  soon  as  new  classes  were  recognised  and 
new  data  are  available,  with  no  need  to  re-creating  the  pre-existing 
modules.  Moreover,  each  module  can  process  data  using  different 
feature  types  and  the  most  adequate  classification  technique  ac¬ 
cording  to  the  different  case  studies  of  natural  and  artificial  rocks. 

Further  work  will  be  devoted  thus  to  the  automatic  classification 
of  other  mineral  inclusions  such  as  feldspars,  micas,  amphiboles 
and  pyroxenes  and  to  different  kind  of  porosity,  extending  the 
modular  classifier  structure.  Actually,  the  use  of  such  an  architec¬ 
ture  may  be  considered  particularly  useful  for  the  heterogeneous 
nature  of  archaeological  materials. 

This  preliminary  work  represents  a  pilot  study  for  the  use  of  the 
proposed  approach  in  ceramic  petrography  and  a  part  of  a  more 
general  methodological  effort  to  classify  pottery  fabrics  by  ANNs. 
Future  work  must  be  then  devoted  to  practical  archaeometric 
application  as  well  as  temper  analysis  and  recognition.  Actually,  a 


Table  5 

Classification  accuracy  of  modular  and  monolithic  classifier. 


Class 

Modular  classifier 

Monolithic  classifier 

Quartz 

Calcareous  aggregates 

Pores 

Total 

Quartz 

Calcareous  aggregates 

Pores 

Total 

Classification  accuracy  (%) 

94.3 

95.7 

79.7 

90.1 

96.1 

90.3 

79.6 

89.0 

Number  of  errors 

76 

45 

235 

356 

52 

106 

236 

394 

272 
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modular  neural  classifier  may  be  used  to  distinguish  temper 
addition  from  natural  occurring  inclusions,  learning  attributes 
characterizing  the  corresponding  petrofacies  or  raw  material 
sources.  At  higher  level,  we  point  to  use  modular  classifiers  to 
distinguish  pottery  fabrics  via  combination  of  compositional  and 
morphometrical  features  of  each  component.  In  this  way,  modular 
architecture  could  be  used  to  attribute  each  sample  to  a  given 
fabric.  For  the  considered  case  study,  this  could  validate  petro¬ 
graphic  observations  which  initially  led  to  distinguish  the  three 
fabrics  as  well  as  sedimentological  interpretations  highlighted  by 
image  analysis  (Eramo  et  al.,  2014). 
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